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Several recent experiments revealed a change of the sign of the first harmonic in the current-phase relation 
of Josephson junctions (JJs) based on novel superconductors, e.g., e?-wave based JJs or JJ with ferromagnetic 
barrier. In this situation the role of the second harmonic can become dominant; in this case, it determines the 
scenario of a 0-7t transition. We discuss different mechanisms of the second harmonic generation and its sign. 
If the second harmonic is negative, the 0-Jt transition becomes continuous and the realization of a so-called (p 
junction is possible. We study the unusual properties of such a novel JJ such as critical currents, magnetic field 
penetration, plasma gap, and microwave response. We also analyze the possible experimental techniques for 
their observation. 

PACS numbers: 74.50.+r, 85.25.Cp 74.20.Rp 

Keywords: Long Josephson junction, sine-Gordon, half-integer flux quantum, semifluxon, 0-pi-junction 



I. INTRODUCTION 

Several novel Josephson junctions (JJs), e.g., 
s-wave/rf-wave JJsL 2 ^, of-wave/if-wave JJs 4 -! 5 -, 
superconductor-ferromagnet-superconductor (SFS) 6 ' 7,8 ' 9,10 or 
superconductor-insulator-ferromagnet-superconductor (SIFS) 
JJ a 11 ' 12 , etcA 3 - can have a current phase relation 

7.v = A t sin((j)), (1) 

with the critical current density j c \ changing its sign, e.g., as 
a function of temperature r 6 ' 7 ' 12 or F-layer thicknes s ^i 11 ! 12 , 
A negative sign of j c \ leads to the formation of a so-called 
K junction with a ground state phase (]) = 7 I 14 ' 15 . Usually, the 
single harmonic current-phase relation (CPR) is adequate for 
the description of the JJ properties; high order harmonic terms 
in j s (§) can be omitted. However at the 0-Jt transition j c \ 
vanishes and the second harmonic becomes important. Then, 
one can write 

7*(<1>) = jd sin((j)) + j c2 sin(2<|>), (2) 

see Ref. and references therein. 

A second harmonic with negative amplitude j C 2 may result 
in the formation of a so-called (p-junction, i.e., a junction for 
which, in the absence of a supercurrent, j s = 0, the phase drop 
across the junction (j) = ±cp such that < cp < Jt, see Ref. [la 

In this paper, we investigate the rich physics of JJs with a 
substantial second harmonic in the CPR, especially when j c 2 
is negative. We discuss short as well as long JJs and propose 
several experimental tests of our predictions. 

The paper is organized as follows. In Sec. [LI] we discuss the 
origin of the second harmonic in the CPR for several types of 
JJs. In the main Sec. [HI] we introduce the model, derive con- 
ditions for the existence of a cp JJs and discuss their properties 
such as critical currents, magnetic field penetration, plasma 
gap, and microwave response. Sec.lIVIconcludes this work. 



II. ORIGIN OF THE SECOND HARMONIC 
A. d-wave superconductors 

The second harmonic in the CPR was observed expri- 
mentally in symmetri c 17 ' 18 45° [001] tilt grain boundary 
(GB) JJs, in tilt-twist GB JJsi2 and in c-axis s-wave/<i-wave 
YBa2Cu307-Nb JJs2°-. The possibility that the first harmonic 
vanishes at some temperature was predicted theoretically 21 
and was observed in experimen t 18 ' 22 . 

Calculations show that the intrinsic second harmonic term 
for bicrystal 45° grain-boundary JJs made of d-wave super- 
conductors is negative 2 —. For more details turn to reviews . 14 ' 23 

Faceting. Even in the cases, when the second harmonic 
was observed experimentally, it was not clear whether it is 
present intrinsically or it is a result of faceting or interface 
roughness. In the latter case, the second harmonic appears in 
the equation for the average (slowly varying) phase, which 
changes on a length scale much larger than the interface 
roughness. The theory of this "effective negative second har- 
monic" was developed by R. Mints and coauthor s 24 ' 25 ' 26 ' 27 
and was also recently discussed in ReffA 6 .. In many cases 
the amplitude of the effective second harmonic \j C 2\ may be 
comparable to or even larger than the amplitude of the first 
harmonic \j c i\- Bicrystal 45° GB JJs made of of-wave su- 
perconductors are one example of such a system where a (p- 
junction may appear. Here, the so-called splintered vortices 
were observed 5 -. 



B. SFS or SIFS junctions 

Next, we consider SFS or SIFS junctions. Often, their CPR 
is sinusoidal only near 7^. At low temperatures the higher 
harmonic terms become more and more important. The cal- 
culations of the CPR in SFS JJs in the clean limit indeed re- 
veal a strongly nonsinusoidal j s (§) dependenc e) 29 ' 30 ' 3 1 . In the 
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dirty limit, in which most experiments are done up to now, 
the js(§) dependence becomes almost sinusoidal when the F- 
layer thickness d F exceeds the decay length of the order pa- 
rameter namely^ 2 - 



d F \ ( d F 



Jd ~ exp 



2d F 



3FI 



(3a) 
(3b) 



In the absence of spin-flip scattering (which is rather unreal- 
istic) in the F-layer of thickness dp, the decay and oscillation 
lengths \p\ and \p% are equal to the characteristic length £, F . 
In the presence of spin-flip scattering 3 - 2 . ^ F1 < <~ F2 , but such 
that ^fi^f2 = Thus, in practice, for dp>%pu the contri- 
bution of the second (and higher) harmonics may be neglected 
everywhere, except for the vicinity of the 0-7C transition. At 
the transition point d F ~ K = §^F2 3> \f\, the first harmonic 
term vanishes and the properties of the SFS junctions are de- 
termined by the second harmonic terrn^. A recent experimen- 
tal study 1 which demonstrates two — > 7C and 71 — > transi- 
tions on the I c (dp) dependence in a SFS J J with a Cuo.52Nio.48 
alloy as an F layer, allows to estimate £, F \ as 1.3 nm, while the 
— > tc crossover thickness d F ' n as 11 nm. This results in a 
very small value of the critical current at the 0-7C transition, 
c.f. Eq. ©, /„ = 7,2 ~ exp(-d$*fei)ln ~ 10- 4 / c0 , where 
I c o ~ exp(— d F ~ K /<^ F i ) is the critical current (dominated by the 
first harmonic) away from the transition, i.e. if there would be 
only the decaying part of the first harmonic without oscilla- 
tions. Therefore, the measured non zero values for the critical 
current at the 0-tc transition (see Refs.^andHH) can hardly be 
explained by an intrinsic second harmonic contribution. The 
intrinsic second harmonic term in SFS JJs at the 0-71 transition 
is positive^ 2 .. 

Inhomogeneous F-layer. Another possible mechanism of 
a negative second harmonic generation in SFS JJs is the inho- 
mogeneity of the F layer thickness near the 0-7C transitions. As 
a result, the JJ consists of alternating and 7C mini-junctions, 
which is similar to the case of faceted GB JJs. The CPR in a 
diffusive SFS junction is described with good accuracy by an 
effective (slowly varying) phase \|/ 



= jd sin(y) +j C 2sin(2v). 



(4) 



For a long JJ (LJJ) with alternating current density (j ca and 
j c b within regions of lengths a and b respectively) one finds 



(a j ca +b j cb ) I '(a + b), 

1 a 2 b 2 (j C a-jcb? 
tfa 24\ j ca \ (a + b) 2 ' 



where 



3>o 



2nfiod' I j c 



(5a) 
(5b) 



(6) 



and /Jod ! as /j()(dj + 2X£) is the inductance (per square) of the 
superconducting electrodes forming the JJ, see Refs. |28| and 



[l6l for more details. The expressions © are valid when a, b <C 
Xj a . When \j c i/2j C 2\ > 1 the ground state of the system has a 
uniform phase (0 or 7t), see below. 

However, in the presence of an applied current a spatial 
modulation of the phase appears with the amplitude 



A(j): 



1 ab(j c 



' Jcb) 



ha 



32|j ca 



sin(\|/) < 1. 



(7) 



This corresponds to the appearance of two types of fractional 
Josephson vortices that are soliton solutions of the double 
sine-Gordon equatio n 34,35 ' 36 . Such vortices were observed by 
SQUID microscopy in 45° grain boundary JJs 5 -. 



C. Experimental techniques 

Let us discuss typical experimental methods and their abil- 
ity to distinguish the sign of the second harmonic. 

1. The most direct technique is a measurement of the 
CPR by embedding the investigated JJ in a SQUID 
loop 17 ' 18 ' 20 . The difficulty here is that not all types of 
junctions can be incorporated into a SQUID with proper 
Pl < 1 values, where = 2I c L/<£>o is the is the induc- 
tance parameter of the SQUID, and I c is the total critical 
current of the J J. 

2. Measurements of I c vs. magnetic field may reveal a 
twice shorter period of oscillations, but the sign of the 
second harmonic is difficult to determine. In fact, our 
simulations show that I C (H) curves for LJJs with strong 
positive or negative second harmonics of the same am- 
plitude look the same for short JJ and qualitatively the 
same for LJJ. 

3. Measurements of sub-harmonic Shapiro steps are also 
quite unreliable to determine the sign of the second har- 
monic, as the situation is very similar to the previous 
case — subharmonic steps appear for both positive and 
negative j c % and look qualitatively the same. More- 
over in JJs with nonvanishing capacitance (Stewart- 
McCumber parameter p c = 2nI c R 2 C/<i>o > 0) fractional 
Shapiro steps appear in any case. Nevertheless, it seems 
possible to extract useful information about the second 
harmonic from the Shapiro step modulation even in the 
presence of a finite capacitance 3 - 7 -. 

4. The presence of a (p-junction is a strong evidence for a 
negative second harmonic. In a LJJ this can be man- 
ifested through the existence of (splintered) fractional 
Josephson vortices of two different kind o 24 ' 25 ' 26,27 that 
are solitons of a double sine-Gordon equation 3 ^. The 
existence or motion of such vortices may be detected 
experimentally using SQUID microscopy 5 . 
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III. RESULTS 

The double sine-Gordon equation which describes the dy- 
namics of the Josephson phase in the LJJ in question is 



tyxx-tytt - [sin(|)+gsin(2(|))] = acj), -y, 



(8) 



where g(T) = ]&/ j c \ is a relative amplitude of the second 
harmonic, which, generally speaking, is a function of the tem- 
perature T. Subscripts x and t denote partial derivatives with 
respect to coordinate and time, respectively. The coordinates 
are normalized to 



4> 



2%fiQ\j c i\d' 



when 



,38 



d' = dj 



( y 



[ y 



(9) 



(10) 



and di ; 2, ^-1,2 are the thicknesses and London penetration 
depths of the superconducting electrodes and d] is the thick- 
ness of the (insulating) barrier. The time is normalized to the 



inverse plasma frequency co pl , where 



C0„i = 



jcl®0 

2nC ' 



(11) 



The parameter a = 1 / ypc is the dimensionless damping pa- 
rameter, and y = j/ j c \ is the normalized applied bias current 
density assumed to be uniform. Further, (j)(x,f) describes ei- 
ther the real phase when the second harmonic is present in- 
trinsically or the phase averaged over facets for the case when 
the second harmonic appears due to faceting. 

The Josephson energy density (per unit of LJJ length) is 
given by 



f/((j)) = £ysg 



n(;'ei)w| 



■f[l 



s(2<|>)]}. 



(12) 



where £/ = 4>o|,/cl|/2Jt sets the characteristic scale of energy 
density, and w is the JJ width. U(§) is defined with accuracy 
of a constant and this constant in Eq. ( fT2l is chosen so that 
U(Q) = 0. This is a natural choice for conventional JJs with 
jcl > and j C 2 = 0, as it corresponds to the energy minimum 
U = which is reached at (j) = 0. In our more general case, 
<)) = does not necessarily correspond to the energy minimum, 
but we still will use the same reference level for t/((j)) to avoid 
confusion. 



A. Ground states and cp-j unction 

Consider a uniform ground state of the system, (j)(x) = 
const. In this case the analysis for the LJJ reduces to the anal- 
ysis of the ground state in a point-like JJ. Let us investigate 
conditions at which one can obtain a cp-j unction 16 , i.e., the 
junction for which 



: for (j)(x) = const. 



(13) 



Substituting expression (O into Eq. ( fT3l we arrive at three 
posible solutions 



(j) = 0; 
<j> = 7t; 
(j) = ±9, 



where 



(14a) 
(14b) 
(14c) 

cp = arccos I — I = arccos I I . (15) 

VdzJ V 2 8j 

The ground state ( 114ct corresponds to the (p-junctioni£. The 
stable solution should correspond to the energy minimum, i.e. 

d 2 U{if) _®odj,($) 



Jcl 



1 



>0. 



(16) 



dd} 2 2k c/(j) 

Substituting each of the solution ( TT~4T > into ( TTol l we obtain the 
stability conditions 

jd > -2j C 2l (17a) 
jd < +2jci; (17b) 



;2 

Jcl ^ n ; 
7TT~ > 2 Jc2- 
2jc2 



(17c) 



Note that in addition to stability condition ( 117c| ) we should 
impose a condition to the argument of the arccos in Eq. ( fTBI ), 



Jd 



2jc2 



< 1. 



(18) 



An analysis of conditions ( I17cb and ( fT8l ) shows that the 9 
ground state may be realized only for 2j C 2 < — \jci\- This 
means that a cp-junction can be obtained from a or a Jt junc- 
tion with a strong negative second harmonic. 

Fig.[TJshows the stable uniform ground states as a function 
of g for jd > (a) and j c \ < (b). The ground state dia- 
gram depends very much on the sign of j c i ■ Fig.QJ shows the 
regions of 0, Jt and ±cp ground states on j c \ - j C 2 plane. 

1. Consider the case j c i > 0, see Fig. [T]i. The state (j) = 
is the ground state of the system for g > —1/2. For g < 
— 1/2, there are two degenerate ground states (j) = ±(p, see 
Eq. <n~3T > . In addition, for g > 1/2 there is a stable state <p = % 
corresponding to a local minimum of energy, but its energy is 
larger than the energy of the ground state (j) = 0. 

2. The case j c i < (7C JJ) is shown in Fig. QJ). The state 
<p = Jt is the ground state of the system for g < 1/2. For 
g > 1/2, there are two degenerate ground states (j) = ±cp, see 
Eq. (TT3T >. In addition, for g < — 1/2 there is a stable state (j) = 
corresponding to a local minimum of energy, but its energy is 
larger than the energy of the ground state <p = 7t. 

One may notice that the ground state diagrams shown in 
Figs. OJ and b are very similar. In fact, one can reduce one 
case to the other by a simple transformation: <p — > 7t — <p, 
g —> —g. In fact, Eq. ([8]l is invariant with respect to this trans- 
formation. In this way, a JJ turns into a Jt JJ, and the bi- 
furcation point turns from g = =Fl/2to±l/2. Applying this 
transformation twice we go back to the initial case. 

Below, without loosing generality, we consider only the 
case jd > 0. The results for the case j c i < can be naturally 
obtained by employing the above mentioned transformation. 
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FIG. 1 : (Color online) Ground state phase as a function of g for (a) 
j c \ > and (b) j c \ < 0, dashed lines indicate stable states with higher 
energy, (c) show the ground state of the system on the plain j c \-j c 2- 



B. Critical current 

When a JJ has two harmonics in the CPR, the natural ques- 
tion which arises is "how the measured value of maximum 
supercurrent (critical current) I c is related to the amplitude 
of both harmonics I c \ and / C 2?" To answer this question we 
rewrite Eq. (O in normalized units 



7(40 = M$)/jei = sin(4>) +gsin(24>). 



Looking for an extremum of y((j)) we find that this extremum 
is reached for the phase (j)o such that 



cos(4>n) = 



-l±v/l+32g 2 



(20) 



Note, that the solution with the "— " sign in Eq. d20b only ap- 
pears for |g| > 1 /2. 
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FIG. 2: (Color online) The dependences y c ±(g) and their asymptotic 
behavior for g — > °°. 

The primary (maximum) critical current is always given by 
substituting the solution ( f20l > with the "+" sign into Eq. (fT9l l 



Jc+(g) 



1 



32|g| 



(21) 

This dependence is shown in Fig. [2] For small g it behaves as 
y c+ (g ^ 0) w 1 + 2g 2 + 0{g% while 



(22) 



For |g| < 1 /2, the critical current is simply given by Eq. ( 1211 ). 
The physics is similar to a JJ with g = 0. 

For |g| > 1/2, the secondary critical current appears, see 
Fig. |2] It corresponds to <|>o with the "— " sign in Eq. ( f20b and 
is given by 



7c- (g) 



1 



32|«| 



1+32^-3 



1 + 32^ + 1 



Note, that y £ _ (g — > ±°°) 
the difference 



(23) 

1 ±g± It is interesting that 



Ajc=Jc+(g)~Jc-(g)~V2, 



(24) 



(19) 



is almost constant as can be seen in Fig. [2] The largest devia- 
tion of 8 % takes place at \g\ = 1 /2. From Eqs. (f2TT > and j23l ) 
one can see that Yc±( — g) = Jc±(+g)> therefore in Fig.|2]we 
show y c ± (g) only for positive g. 

An important question is: can one observe y c _ in experi- 
ment? To answer this question we have to consider the case 
of negative and positive second harmonic separately. 
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FIG. 3: Josephson energy Uj of the system as a function of the phase (]) for three cases: (a) g = — 1, (b) g = and (c) g = 1 



1. For g > 1/2, and at y = the system has two stable 
states (j) = (with lower energy) and (j) = K (with higher en- 
ergy), see Fig. 0c). The critical current y_ corresponds to 
the "depinning" of the Josephson phase from the high energy 
71-state, while y c + corresponds to a depinning from the low 
energy 0-state. If initially the junction is in the 0-state, one 
just measures y c+ (1211 1. On the other hand, if initially the JJ 
was in the 7i-state, at y c _ the Josephson phase starts moving. 
Depending on damping, it may either result in a stationary 
phase motion (for low damping) or the phase may go down 
to the neighboring lower energy 0-state and stick there (for 
high damping). Thus, an underdamped JJ which was initially 
in the Jt-state will switch to the resistive state at y c _, while 
the overdamped JJ will just switch from the 7i-state to the 0- 
state at y e _ and will switch to the resistive state only when the 
bias current is further increased above y .+. We note here that 
the probability to find the JJ initially in the 7t-state may not 
be very low. For example, for large g, the ratio of the energy 
difference between the 0- and the K states to the maximum 
barrier height goes like [U(n) — U(Q)]/U max « l/g, i.e., be- 
comes negligible. This means that during switching from the 
voltage state to the Meissner state the phase may stick in the 
7i-state with a probability close to 50%. If one is able to deter- 
mine experimentally both / c _ = / c iy c - an d h+ = Icijc+, one 
will then be able to calculate I c \, I c % and g from experimental 
data. For large g this calculation is straightforward: 



Ic2 



^(Wc-); 
hi 



(25) 
(26) 



2. For g < —1/2, the system has two degenerate stable 
states (j) = ±(p = ±arccos(— l/2g) (cp-junction), as shown in 
Fig. |5J a). The critical current y ._ corresponds to the escape 
of the phase from the state — (p towards the state +(p, while 
y c + corresponds to the escape of the phase from the state +(p 
towards 27T — cp over the large potential barrier. If initially the 
JJ is in the +(p state, in experiment one observes only y + . If 
initially the system is in the — (p state, then, upon exceeding 
y_, the phase moves towards the +(p-state and either ends 
up being trapped there (typical for an overdamped JJ) or may 
continue moving further switching the JJ into the voltage state. 
Which of these possibilities is realized depends not only on 



damping but also on the height of the potential barrier which 
the phase should overcome to keep moving continuously. This 
barrier, in turn, depends on g, see Fig. |4] 
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FIG. 4: The Josephson potential energy U (0) at y = y e _ (g) for g < 
go g = gc~ -1-26 and g > g c . 

A numerical study shows that the JJ with vanishing damp- 
ing will not be trapped in the +tp-state for g < g c ■ « — 1 .26. At 
this value of g the height of the barrier is equal to or smaller 
than the initial energy of the phase-particle, as shown in Fig. [4] 

Thus, for a sufficiently large and negative second harmonic, 
one can observe two critical currents on the I-V characteristic 
(IVC). Actually, the typical way to trace the IVC by sweeping 
the current / from — 7 max to +/ max an d then back to — 7 max will 
most probably result only in ±/ c _ to be visible. To see +I C + 
one should sweep from to +/ m ax, to 0, to +I max , etc., i.e. 
only the positive part of the IVC. To see —I c + one should trace 
only the negative half of the IVC. Simulations of the IVCs 
using different sweep sequences confirmed this prediction. 



C. Plasma waves 



One of the unique properties of Josephson junctions is the 
presence of the plasma frequency <s> p . In a short JJ, co p is the 
frequency of the eigenoscillations of the Josephson phase at 
zero applied bias current. In a LJJ this is true only if the phase 
is uniform. In general, electromagnetic waves can propagate 
along the LJJ only if their frequency is above co ;) . The disper- 
sion relation (£>(k) = \/\ +k 2 in a conventional LJJ has a gap 
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from to 1. In quantum circuits the plasma gap h& p ^> IcbT 
protects the circuit from thermally excited plasmons. It also 
defines the attempt frequency during the thermal escape from 
the zero voltage state as well as it defines the energy level 
spacing in the quantum regime. Therefore it is important to 
look into the dispersion relation and the plasma gap in the LJJ 
with a second harmonic in the CPR. 

The derivation of the dispersion relation co(£) for plasma 
waves propagating along the junction is straightforward. We 
substitute the small amplitude wave solution (j) = § s (y) + 
Aexp[/(A:x — cof)] into Eq. ([8]) without a damping term. <j) s (y) 
is the uniform phase, which is a solution of the static Eq. (|8). 
When y = 0, it takes the value of (j). s = J14at for \g\ < \, 
§ s = 0,71 ( I14all4bb for g > \ and = ±ip ( I14cb for g < -\. 
As a result one gets the dispersion relation 

>(*) = ^ p (g,y)+k 2 , (27) 



CO 



where (O p (g,y) is the plasma gap, which depends on the 
ground state, the amplitude of the second harmonic g and on 
the applied bias current y. In a conventional (g = 0) LJJ 



<(0,Y) 



(28) 



where the superscript stands for the ground state (j) = 0. In 
the general case S) p (g,y) can be calculated as 



<O p (g,y) = a/cos [<k(y)]+2gcos [2(j). $ (y)]. 



(29) 



Since the static phase § s (y) is a solution of the transcendential 
equation sin^ + gsin2(]) s = y, one cannot write the explicit 
expressions for & p (g,y) except for some limiting cases. When 
y = one obtains 



y/2g+l, in 0-state, g > 



1 



\/2g— 1 , in 7C-state, g > +-; 



/ 1 1 
— -2g, incp-state, g< --. 



(30a) 
(30b) 

(30c) 



The plot of the gap (£> p (g,Q) is shown in Fig. [5] The plasma 
gap closes and opens again at g = —1/2. The slope from the 
left is — 2e, while from the right it changes abruptly to +4e, 
where e= \g+ j|. 

The dependence of the plasma gap (O p (g, y) on y looks qual- 
itatively similar to Eq. d28l ), but with properly renormalized 
(Q p (g,0), see Eqs. (f30|, and y c (g), see Eqs. (ED and d23l 

For the weakly biased y <C 1 state corresponding to the 
ground state (j) = one finds 



<(g,Y)«V2l+l- 



1 



IT 2 , 



0. 



4(2^+l)V4 

For y — > y c+ , such that 8y = y c+ - y <C 1 the result is 

1/8 



(q + 3)(q-l)q 2 



16g 2 



(31) 



(32) 
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FIG. 5: The gap in the linear plasma wave spectrum as a function of 



where we have introduced the new variable 



<?= v / l+32g 2 , 



(33) 



which is also used below to make formulas more compact. 

For the ground state = Jt, the weakly biased state will 
have the plasma frequency 



«FJg,y)^^2g—\- 



-1 



FT 2 , 



4(2* -1)5/2' 
while in the pre-critical region 5y = y _ (g) — y <C 1 

1/8 



(34) 



<(S,Y) 



-3)(*+iy 



(35) 




-0.5 0.0 0.5 1.0 1.5 2.0 
bias current, y 

FIG. 6: Op{g,y) as a function of yfor g = — 1 and the ground state 
(J) = +(p (thick solid line). Other lines show the approximation of this 
dependence in the vicinity of — y c _ (dash dotted line), y = (dashed 
line) 

A much more interesting behavior of tt>p(g,y) is observed 
for (p-junctions. The behavior of C0p(g, y) as a function of y is 
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shown in Fig. [6] One can see that (af,(g,y) has a maximum not 
at y = 0, but at 



(36) 



2 1? 



2,? 



2Vl-4/ y 8V(l-4g 2 ) 5 



(32s 2 4 
(37) 

For y in the vicinity of — y e _ , we define 8y = y + y ■_ and 
obtain 



(?-3)(?+l)9 



21 



16g 2 



1/8 



(38) 



For y in the vicinity of y + , we difine 8y = y £ + — y and obtain 



(q + 3)(q-l)q 2 



1/8 



'Sy. 



(39) 



D. Josephson Vortices 



The Josephson vortices which may exist in LJJs with a sec- 
ond harmonic can be very different from Josephson vortices 
in conventional LJJs. For the case \g\ < A they still resemble 
usual sine-Gordon kinks. For g < — i there may be two kinds 
of solitons in the system: the smaller one (bj (x)carrying the 
flux <t>i = <t>Q<p/n (2cp vortex) and a bigger one fai*) carry- 
ing the flux 4>2 = *o - *i = <3>o(l — <p/t) (2;t — 2<p vortex), 
see Refs i 16 ' 24 ' 25 ' 26 ' 27 for details. The shape of these vortices is 
given by the following formulas 



<t>iM 



2 arctan 



, xsincp / (p 
tanh I — — - tan ' 



§2(x) = 7C + 2 arctan 



tanh[ ^cot(- 



(40a) 
(40b) 



, <j> (b 
: 4 sin — cos — . 



(41) 



From here, by multiplying with (j) v and integrating we obtain 

2 



d§ 

dx 



-2cos(b — sin <b + C, 



(42) 



for the ±(p state, accordingly. 

This feature is a direct consequence of the asymmetry of 
the potential well in which the phase-particle is trapped in the 
ground state. For comparison, this well is symmetric for and 
71 ground states. 

The behaviour <x>p(g,y) at y— > can be approximated as 



where C is an integration constant. To satisfy the boundary 
conditions at x — > ±°° (vanishing magnetic field (b v ) one must 
set C = 2, so that the r.h.s. of Eq. (l42b becomes equal to 
4sin j, i.e., 



Integrating we get 
13)^. 



and is shown in Fig. [7] 

Such kinks were extensively studied in the framework of 
the double sine-Gordon equation and have some very inter- 
esting properties, e.g., one of them has an eigen oscillation 
frequenc y 35 ' 39 . 

Here, we will focus our attention on the point g = — | where 
the transition from to (p-junction takes place. Let us write a 
static version of Eq. (O, i.e., the analog of a Ferrell-Prange 
equation (g = — ,) 



(b(jt) = ±2 arctan t— r +2%n, , 
\x\ 



(43) 



(44) 



where the ± sign corresponds to solitons of positive and 
negative polarity. Note that if one takes the + sign and n = 
for x < 0, then one should take the — sign and n = 1 for 
x > in Eq. (PPfl i. One can see that the soliton tail is non- 
exponential, namely <b w l/x for x — > ±°°. 



E. Penetration of magnetic field 

Let us consider the penetration of magnetic field into a 
semi-infinite LJJ (x — 0. . .°°). Note that Xji (0 just defines 
the unit length. It does not have the sense of the magnetic field 
penetration depth anymore (see below). 

The phase distribution over the junction is determined by 
the static version of Eq. © 



tyxx = sin(b + ^sin(2(b) 



(45) 



Below we consider the penetration of magnetic field for 
three different ground states of the system. 

a. The ground state (b(°°) = 0, which is realized for g > 
— 1/2. In this case the first integral of Eq. d45l l is 



i<b 2 = (l-cos^) + |(l-cos2^). 

The magnetic field penetration depth can be defined as 
1 



(46) 



H(0) Jo 



H{x) * = ^w-m] =hi [*H-m] 

v ' 2%AH(0) 



h(0) 
(47) 

where the dimensionless field is defined as usual, h(x) = 
§ x (x) = //(x)27cAA,yi/4>o, and 

A = <// + ki tanh +)i2tanh^J^ , (48) 

is the effective magnetic thickness. 

The boundary condition is h(Q) = h. From Eq. ( f46b we get 



■4(2g+l)y + h 2 =0, 



(49) 
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FIG. 7: The shape of the small and the large fractional vortices (140b . The Josephson energy profile U (<|)) is shown in the right plot and also as 
a grayscale background in the left plot. 



where we defined y = sin 2 [0(0) /2]. Solving for y we get 



Substituting (j)(0) 



(2g+l)- v /(2g + l)2-2g/z2 

2arcsin(,/y) into Eq. ( f4Tb , we find 



A-r, — 



/i< l, 



(50) 



(51) 



where the superscript stands for the ground state (j) = 0. At 
g = — A the screening length diverges (screening is non expo- 
nential) and X° H =Xj\ */| for our case of small field A<1. 

Let us also find the penetration field h p , i.e. the field at 
which vortices start penetrating into the LJJ. From Eq. d49b 



h 2 (y) = -8gy 2 +4(2g+l)y. 



(52) 



We should find the maximum value of h under the constraints 
that < y < 1. For |g| < 1/2, the maximum value h° p = 2 
is reached for y = l. For g > 1/2 the maximum value h® — 
(2g + 1 )/y/2g is reached for y=(2g + l)/4g. 

b. The state 0(°°) = K, which is realized for g > 1/2. In 
this case the first integral of Eq. d45l l is 



i0 2 = (-l-cos0) + |(l-cos20). 
Following the same procedure we arrive at 



A< 1. 



(53) 



(54) 



Note, that A-^ diverges when g — ► 1/2 + 0. 

To find the penetration field we have to search for the max- 
imum of 



h 2 (y)=-Sgy Z + 4(2g+l)y-4, 



(55) 



with the constraints that < y < 1. For g > 1/2 the max- 
imum is reached at y = (2g + 1 ) /4g and is equal to hJi = 
(2g-\)j^Tg. 

c. The ground state (j) = ±cp, which corresponds to the 
domain g < —1/2, The first integral reads 



1 



as 



(coscp — COS0)^ 



(56) 



2 " 2cos<p 

At x = 0, A - = h and following the same procedure as before 
we get 



<p V C0S( P 

^-h = Ml — : 

sincp 

Forg^ -1/2-0, Xl->\n 



= X 



-2g 



(57) 



16 



V-2(2S+1) 



in accordance with 



the previous result 

To calculate the penetration field hp we note that the LJJ 
may be in one of the ground states ±(p. Without loosing gen- 
erality we assume that it is +(p. We also assume that the LJJ is 
very long, but still has two edges at x = ±L/2. From Eq. d56l ), 
at the left edge of the LJJ (x = —L/2) we have 



(coscj) — coscp) 



(58) 



^/coscp 

and the maximum field (penetration field) h^ pl corresponds to 
the penetration of the 2(p-soliton (I40at into the LJJ from the 
left edge. 

^ = ^^ = -^, (59) 
p [ y/cosy V^2g 

Similarly, at the right edge x = +L/2 we have 

(coscp — coscb) 
& = k=- ^== , (60) 



^/cos9 
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and the maximum field (penetration field) h% corresponds to 
the penetration of the (2k — 2(p) -soli ton ( I40bl i into the LJJ 
from the right edge. 



and cos (02)^ 



<p _ (coscp 

h p2 - 



1) 



1-2* 



^/cosf ^f=2g 



(61) 



Since < h^ )2 , the penetration of the flux in terms of vor- 
tices into the junction is not symmetric. At h < hJ A there are 
no vortices inside the LJJ, but the tails of two different kinks 
penetrate into the LJJ from the edges: a tail of a 2(p-kink from 
the left edge and the tail of a (2k — 2(p)-kink from the right 
edge. At h ^ < h < h^ 2 , a single 2(p-kink enters the LJJ from 
the left side and the phase behind the kink (between the kink 
and the left edge) sets to — (p. Thus, the next kink that should 
enter into the junction from the left side is a large (2k — 2(p)- 
kink. But the field h is still below ht 2 — the penetration field 

for the big kink. Thus, upon exceeding hj n , a single small 
vortex enters the LJJ and the field does not penetrate further 
until the field exceeds h 2 . Upon exceeding h^ 2 , a chain of al- 
ternating big, small, big, small, etc vortices enter into the LJJ 
from the left and from the right ends until they fill the whole 
LJJ with a certain density. 

Fig- E visually summarizes the results on A.//(g) and h p (g) 
obtained above. 



F. Dependence I C (H) 

The calculation of the I C (H) dependence for a short JJ of 
the length 2a = L and the width w follows the usual procedure. 
We assume that the phase is linear with magnetic field, 



= fcc + 0o, 



(62) 



where /; = 2H/H C \ is the dimensionless applied magnetic 
field. It is normalized to H c \ = Qq/kKjiA which is the pene- 
tration field of the JJ with g = 0. 

The total supercurrent through the junction normalized to 
the maximum supercurrent of the first harmonic j c \2awk?j is 
given by 

1 f +a 

i s (h,g,§ ) = — / [sin(/ix + o ) +§sin(2/jx + 20o)] dx 
2a J - a 

2 sin (ha ) sin (0o ) + g sin (2ha ) sin ( 20o ) 
2ha 



cos(0i) = 



COS (02 ) 



-1 + y/\ +32g 2 cos 2 (ha) 

Sgcos(ha) 
1 - y/\ + 32g 2 cos 2 (ha) 

8gcos(/ia) 



(65) 



(66) 



To get I C (H) we substitute 0i.2 from Eqs. (l65T l and (l66l l into 
Eq. (l63l . i.e., I c {h,g) = /(/z,g,0 1,2 (/?,£))• in principle, we get 
two branches of I c (h). 

The solution given by Eq. d65l l always has |cos(0i)| < 1 
and, upon substitution into Eq. 1631 . yields the upper branch 
I C+ (H), which corresponds to I C (H) in a conventional JJ with 
g = 0. Instead, the second solution given by Eq. (|66l l has 
I cos(02) I > 1 if \g\ < 1/2, i.e. it is irrelevant for \g\ < 1 /2. For 
\g\ > 1/2, |cos(02)| may or may not be below one, depend- 
ing on the value of the applied field, and, for some ranges of 
magnetic field, yields 7 C _ (H). The calculated I c ±(h) plots are 
shown in Fig.|9]for different values of g. One can see that fea- 
tures with twice shorter period in H appear as \g\ increases. In 
the limit \g\ — > °°, I C+ (H) again becomes a Fraunhofer pattern, 
but with twice shorter period. 

We note that I c (h,g) looks exactly the same as I c (h,—g). 
Indeed, changing +g to — g changes the sign of cos (0 12) in 
Eqs. ((65) and (1661 . This corresponds to the change from 
0i ; 2 to K — 01,2- Looking at Eq. d63l one sees that sin(0o) 
is a symmetric function of 0o with respect to 0o = k/2, i.e., 
does not change, while gsin(20i.2) = — gsin(2(7C — 0i,2)) i.e. 
I s (h,g,tyi l2 ) =/ s (/i,-*,7C-0i, 2 ). Thus, I c (h,g) =I c (h,-g). 
From a practical point of view, this means that one cannot 
distinguish between a J J with positive and negative second 
harmonic by measuring I C (H) in a small JJ. 

For \g\ < 1/2, see Fig. |9|a), the plots I c +(h) look very 
similar to the usual Fraunhofer dependence except for the 
region around h = 0. In practice, one can detect the pres- 
ence of a weak second harmonic by measuring I C (H) and 
calculating the ratio of I C (H\) at the first side maximum to 
the / c (0). In a conventional JJ this ratio is w 2/371 « 0.21. 
In the presence of the second harmonic I C (H\) stays almost 
the same, while / c (0) changes with g, see Eq. (fJTJ, reaching 
3^/3/4 w 1.3 at \g\ -> 1/2. Thus l c (Hi)/l c (0) changes from 
0.21 to8\/3/277Cw0.16for \g\ from to 1/2. 

For |g| > 1/2, see Fig. |9jb)-(c), the second lower branch 
I c -(h) can be seen. Note that this branch appears not for all 
values of h, but only for those where the r.h.s. of Eq. ( |66l l < 1 . 
The values of the r.h.s. of Eq. d66l > are shown in Fig. |9jb)— (c) 
by the dashed line. Note, that 7 C _ (h) < I c +(h) for all values of 
h. 



Note, that ha = K<t>/<pQ is a normalized flux inside the junc- 
tion. To find the maximum value of the supercurrent we 
should find the maximum of 7 s (/z,0o) with respect to 0o. The 
maximum is reached when 3/ ( /30o = 0. This yields the fol- 
lowing equation for cos(0o): 

8^sin(2/ia)cos 2 (0o) +4sin(/ja)cos(0o) — 4gsin(2/ia) =0. 

(64) 

This equation has two solutions, which we denote as cos(0i) 



G. Shapiro steps 

A junction with a second harmonic shows not only integer, 
but also semi-integer Shapiro steps when a microwave voltage 
is applied, even when the capacitance of the JJ C = 0. Here 
we consider a simple model of applied dc + ac voltage 



V(t) = V + Vicos(a>f). 



(67) 
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FIG. 8: Dependence of (a) the small magnetic field penetration depth 'kfj(g) and (b) the vortex penetration field h p (g) on g. 



From the second Josephson relation 0, = 2kV /<t>o, the phase 
in the junction will have the form 

0(0 = fya + Gjt + Asin(coO, (68) 

where 0o is an arbitrary phase shift (integration constant), 
CO/ = 27rt / o/ < t > o is the Josephson frequency (frequency of in- 
ternal oscillations in the JJ), and A = V\ /co is the amplitude of 
the phase oscillations. 

From here on we follow the standard procedure 4 -^ We 
substitute the phase ansatz d68l ) into the expression for 
the supercurrent ( fT9l l, expand terms like cos(Asin(cof)) and 
sin(Asin(cof)) in terms of Bessel functions, and find the dc 
component e.g. by averaging over time. Below we present 
the results for CO/ = nco, where n is an integer or half-integer 
number. 

For integer n the general formula for the supercurrent 
reads 4 -^ 



It is obvious that the maximum contribution to the dc current 
takes place for sin(20o) = ±1, therefore-' 

i^{A)=g\J 2n {2A)\. (72) 

As we see, the semi-integer Shapiro steps appear as soon as 
g t^O. The behavior of yJJ 13 " (A) for« = |is shown in Fig. [10] 
for \g\ = 1, but can be scaled to get the dependence for arbi- 
trary \g\. 

Recently, a formula for integer and semi-integer Shapiro 
steps which takes into account the capacitance of a JJ was 
obtained 37 . This formula is valid in the high frequency limit 
and, in principle, it allows to determine the sign of the sec- 
ond harmonic experimentally. On the other hand, this formula 
does not take into account the possibility of chaotic dynamics 
in the underdamped system, so its application may be prob- 
lematic in some cases, and one should rely on numerical si- 
multions. 



Y(A,0 O )„ = |sin(0 o )/«(A) +gsin(20 o )/2„(2A)| . (69) 



H. Zero Field Steps 



As we see the contribution of the Josephson supercurrent into 
the total dc current depends on 0o. The amplitude of the 
Shapiro step corresponds to the maximum of y(0)„ with re- 
spect to 0o- The maximum is reached for 



COS(0 O ) 



J„(A)±^(A) + 32g 2 jj,(A) 

8ghn(A) 



(70) 



Each of these two solutions is relevant only if |cos(0o)| < 1. 
If, for given A, |cos(0o)| < 1, we calculate 0o as arccos of 
Eq. (|70T i and substitute it into Eq. ([69). We get a rather bulky 
expression for yJf ax (A), which we do not show here. Note, 
that taking — arccos () of Eq. (f70l produces the same result. 
We also note that , yJJ lax (A) does not change if we change +g to 



For semi-integer n, the supercurrent is given by 41 



y(A,0o)„=§sin(20o)7 2 „(2A). 



(71) 



The presence of mobile fractional vortices for g < — 1/2 al- 
lows the observation of half-integer zero field steps (ZFS) on 
the current-voltage characteristic similar to classical integer 
ZF S 40 ' 42,43 ' 44 . The dynamics can be described as follows. A 
small vortex situated inside the junction subject to a driving 
force due to the bias current, moves along the junction. When 
it arrives at the edge x = L, the boundary condition requires 
that 0. X (L) = 0. As in the analysis of integer ZFS, we satisfy 
this boundary condition by considering a collision of a small 
2tp-vortex with its "image" — a small — 2cp-antivortex — sit- 
uated outside the LJJ. Thus, instead of treating the collision 
of a vortex with the boundary at x = L, we treat the collision 
of a vortex with an antivortex at x = L in the absence of the 
boundary. Such a collision of a small vortex and an antivortex 
should inevitably result in the appearance of a large (2k — 2(p)- 
vortex and a (2cp — 27t)-antivortex after collision, see Fig.fTTI 
Thus, a 2tp- vortex colliding with the boundary, reflects back as 
a (2cp — 27i)-antivortex. Then the bias current pulls this large 
(2cp — 2k) -antivortex towards the opposite boundary at x = 0. 
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FIG. 9: (Color online) The dependences I c + (ti) and / c _ (/i) calculated 
using Eq. j 63 1 > with (j)] 2 given by Eqs. ((65) and l |66l >. (a) for g = 0, 
±0.3, ±0.5 only the main branch I c +(h) is present; (b) for g = ±1; 
(c) for g = ±3. In (b) and (c) the r.h.s. of \66\ are shown by the 
dotted line. The solution corresponding to I c -(ti) exists only when 
this line goes below 1 . 



There, it reflects as a small 2(p-vortex, which is pulled towards 
the boundary x = L, and so on. The flux transfer per period 
is <&i — (— <t>2) = Oo, therefore the voltage across the junc- 
tion is V = <t>Qu/2L. When the bias is icreased the velocity of 
vortices u — > co and the voltage approaches V — > Qqcq/2L i.e. 
half of the usual ZFS. 

However there is a problem. Even in the lossless case, the 
reflection of a small vortex as a large vortex requires an addi- 
tional input of energy (which can be taken from kinetic energy 




2 4 6 8 10 12 14 16 18 20 
microwave amplitudes [norm, units] 



FIG. 10: (Color online) The behavior of y^' dx (A) for n = \, §, f . 



or from the current source) to compensate for the energy dif- 
ference between the small and the large vortices. Therefore, 
such a reflection may take place only when <t>i « <E>2, i.e. for 
rather large \g\. Indeed numerical simulations show that half 
integer ZFS can be observed for g < — 1. The range of bias 
currents in which the first half-integer ZFS is stable is shown 
in Fig.Q~2]for different values of g and a = 0.1. On one hand, 
the value of a should be smaller than 1 to observe the dynam- 
ics. On the other hand, if a is too small, upon reflection one 
excites too many plasma waves that destroy the solitons and 
switch the junction into the normal state. 

Note that the usual integer ZFS exist in such a LJJ also, 
and they are much more stable. The shape of the soliton is of 
course different from the classical sine-Gordon kink (having 
two maxima in <j> x (jc), but its general behavior is the same. 



IV. CONCLUSIONS 

We have investigated some properties of long Josephson 
junctions with an arbitrarily strong amplitude y' C 2 of the sec- 
ond harmonic in the current phase relation, which is always 
the case in the vicinity of the 0-71 transition where the ampli- 
tude of the first harmonic j c \ vanishes. We have shown that in 
the case of low damping and \g\ > 1/2 (g = ja/jd) one can 
observe two critical currents and determine the amplitude of 
the first and the second harmonics in the current-phase rela- 
tion experimentally. 

For \g\ > 1/2 the second critical current may also be seen 
on the I C {H) dependence, while the dependence itself devel- 
ops twice shorter periodicity in H. Unfortunately, one can- 
not deduce the sign of the second harmonic from the shape 
of the I C (H) dependence since it depends only on \g\. For 
\g\ < 1/2 the I C (H) dependence changes only slightly, but the 
presence of the second harmonic can be noticed by compar- 
ing the heights of the principal maximum and the first side 
maxima. 

In the presence of the second harmonic, half integer Shapiro 
steps appear. Their amplitude is given by the simple formula 
(|72~t . The amplitude of the integer Shapiro steps depends on 
the rf amplitude in a complicated way. Unfortunately, one 
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FIG. 11: Reflection of a small fractional vortex <£>] from the edge of the LJJ as a large fractional antivortex — <1>2- The junction is situated at 
x < 0, the edge is at x = 0. The reflection can be represented as a collision of a small vortex +4>j with a small antivortex — <I>i at point x = 
in an infinite LJJ. They turn into a large antivortex — <E>2 and a large vortex +$2 after collision. Curves 1 to 4 show the evolution of the phase 
with time during this process. The Josephson energy profile U ((])) is shown in the right plot and also as a grayscale background in the left plot. 



cannot deduce the sign of the second harmonic from the mod- 
ulation of the Shapiro steps amplitude. 

The most interesting property of the investigated system is 
that for g < —1/2 two types of fractional vortices may exist 
in the LJJ. One of them carries the flux <t>i and the other <t>2 
such that <t>i + <t>2 = 't'o ■ The smaller vortex (carrying the flux 
<t>i < <t>2) penetrates into the LJJ easier than the larger one, 
resulting in two different penetration fields H p \ and H p 2- The 
presence of the fractional vortices and, therefore, of a strong 
negative second harmonic, may be detected experimentally by 
the observation of half-integer zero field steps that should ap- 
pear in LJJs with moderate damping 0.01 < a < 1. 
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